Effect of an impulsive force on vortices in a rotating Bose-Einstein condensate 
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The effects of a sudden increase and decrease of the interatomic interaction and harmonic-oscillator 
trapping potential on vortices in a quasi two-dimensional rotating Bose-Einstein condensate are in- 
vestigated using the mean-field Gross-Pitaevskii equation. Upon increasing the strength of interac- 
tion suddenly the condensate enters a nonstationary oscillating phase which starts to develop more 
vortices. The opposite happens if the strength is reduced suddenly. Eventually, the number of vor- 
tices attains a final value at large times. Similarly, the number of vortices increases (decreases) upon 
a sudden reduction (augmentation) in the trapping potential. We also study the decay of vortices 
when the rotation of the condensate is suddenly stopped. Upon a free expansion of a rotating BEG 
with vortices the radius of the vortex core increases more rapidly than the radius of the condensate. 
This makes the counting and detection of multiple vortex easier after a free expansion. 

PACS numbers: 03.75.Fi 
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I. INTRODUCTION 

Since the successful detection of Bose-Einstein con- 
densates (BEC) in dilute trapped bosonic atoms at ultra- 
low temperature, one problem of extreme interest is the 
formation of vortices in a rotating condensate in an 
axially symmetric trap. The experimental observation 
[H, H ^, H, H, H, ^ of vortices in the condensate 
has intensified theoretical investigations on vortices us- 
ing the mean- field Gross-Pitaevskii (GP) equation 
For large frequencies of rotation, a very large number 
of vortices arranged in a triangular lattice, frequently re- 
ferred to as a vortex lattice, have been observed. The vor- 
tex lattices are highly excited collective states of BECs 
with an angular momentum of up to 60?i per particle 
1^. There have been many theoretical studies on dif- 
ferent aspects of BEC [l||l| and specially, on the 
vortices in a rotating BEC in axially symmetric traps 

ill, 
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|2^, |2^, |24|, |2^, |2e 

When an ordinary fluid is placed in a rotating con- 
tainer the fluid rotates with the container because of vis- 
cous force. For superfluid ^He(II) in a rotating container 
no motion of the fluid is observed below a critical rota- 
tional frequency. Above the critical frequency, lines of 
singularity appear in its velocity field. These lines of sin- 
gularity correspond to vortices and appear in the form of 
a lattice specially for large rotational frequencies. Super- 
fluidity exhibits many peculiar features: absence of vis- 
cosity, reduction of moment of inertia, occurrence of per- 
sistent current, appearance of quantized vortices among 
others |l3[ . The formation of vortex lattice in a rotat- 
ing '*He(Tl) is intimately related to its superfluidity and 
is a manifestation of quantized vortices ^, ^ ^, ||, ^ . 
However, because of the strong interaction between ^He 
atoms, the theoretical description of this system is not 
an easy task. 

Similar vortices can be generated in theoretical mean- 
field models of trapped BEC based on the GP equations 
[15[ 1^, 1^, |l[ Mi^nd have 

been observed experimentally |E[, ffB,^ pTRiWH] in BEC 



with repulsive interaction. In contrast to liquid ■*IIe(II), 
a trapped BEC is a very dilute and weakly interacting 
system, which makes a mean-field analysis appropriate. 
This analogy of a trapped BEC with liquid ^He(II) sug- 
gests the presence of superfluidity in the BEC. Many dif- 
ferent ways for generating quantized vortices in a BEC 
have been suggested , e.g., via spontaneous vortex 
formation in evaporative cooling p9[ |, via controlled ex- 
citation to an excited state of atoms by stirring a 
BEC using an external laser with an angular frequency 
above a critical value f or g enerating a quantized vortex 
with a single quantum |15| , or by the rotation of an ax- 
ially symmetric trap with an angular frequency above a 
similar critical value |22[ . 

Experimentally, quantized vortices have been gener- 
ated by different groups in a rotating BEC using a vari- 
ety of techniques. Vortices have been detected in a *^Rb 
condensate in a cylindrical trap by Madison et al at ENS 
where angular momentum is generated by a stirring 
laser beam. At MIT Raman et al Q also studied the nu- 
cleation of vortices in a Bose-Einstein condensate stirred 
by a laser beam. They have been observed by Matthews 
et al at JILA in coupled BEC's comprised of two spin 
states of ®^Rb in a spherical trap, where angular momen- 
tum is generated by a controlled excitation of the atoms 
between the two states. The tilting of the vortex axis in a 
nonspherical trapping potential has been studied at JILA 
by Haljan et al [||. Vortices have also been generated by 
evaporatively spinning up a normal gas of ®^Rb and then 
cooling below quantum degeneracy at JILA by Haljan et 
al [0 . A large number of vortices in the form of a lattice 
have been observed at MIT by Abo-Shaeer et al 0, ||] 
as well as at ENS by Madison et al Q for sufficiently 
large rotational frequency. The interesting dynamics of 
the formation and decay of vortices has also been stud- 
ied 1^. Madison et al |2) also confirmed that before a 
vortex pattern is formed the condensate passes through 
dynamically unstable distorted configurations. The dy- 
namical instability transforms the distorted condensate 
to an appropriately symmetric condensate with vortices 
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There have been theoretical studies of how vortices are 
generated in a BEC via these dynamically unstable con- 
figurations using the mean-field GP equations in two | |2^ 
and three |24 dimensions. The theoretical study of the 
dynamical instability using mean-field approach in the 
generation of vortex is of interest as one can then com- 
pare the results of the models with experiments. This will 
provide a more stringent test of the mean-field theories 
than the comparison of stable configurations of BEC. 

However, there is another type of instability which has 
been studied in relation to a nonrotating BEC under the 
application of an impulsive force. One can now suddenly 
change the trapping potential on a stable BEC. This can 
be done in experiment by varying the current in the coil 
responsible for confinement. One can also change the 
strength of interatomic interaction suddenly b y c hanging 
the magnetic field near a Feshbach resonance |g8[ . These 
set in dynamical instabilities in the BEC, which can be 
studied both experimentally |^ and theoretically |30t] . 

Similar impulsive force can be applied on a rotating 
BEC with vortices and the resultant oscillation studied. 
One can study how the number of vortices change after 
the application of the impulsive force on a rotating BEC. 
In this paper we perform such mean-field analyses of the 
dynamics of a rotating BEC subject to three types of 
impulsive forces. We consider a sudden change in (i) the 
strength of interaction near a Feshbach resonance 
induced by a variation of an external magnetic field, (ii) 
the trapping potential, and (iii) the rotation of the BEC. 
Specifically, we study the dynamics when the strength 
of interaction or the trapping potential is suddenly in- 
creased or reduced. We also study the decay of vortices 
once the rotation of the condensate is suddenly stopped. 
There have been certain experimental measurements on 
such decay of vortices [|[ |[. Finally, the trapping po- 
tential can be removed suddenly in a rotating BEC with 
vortices and the resultant expansion studied. 

We base our study on the numerical solution |l^, 
p3| of the nonlinear time-dependent mean-field GP equa- 
tion jl^, which should provide a faithful description of 
the formation and evolution of vortices. The experiments 
on the vortices used a highly anisotropic trap which may 
be simulated by a quasi two-dimensional model. Hence 
as in the recent studies on the subject in Refs. 
we employ the GP equation in two dimensions. We use 
a split-step method where in the time evolution of the 
GP equation the kinetic energy and angular momentum 
(derivative) terms in the two directions are dealt with in 
two independent steps by the Crank-Nicholson rule com- 
plimented by the known boundary conditions |p3| , 
All the harmonic oscillator and the nonlinear terms are 
treated in a third separate step. We find that this time- 
dependent split-step approach leads to good convergence. 

In Sec. II we describe briefly the time-dependent GP 
equation for a rotating BEC and a numerical method for 
its solution. In Sec. Ill we report the numerical results of 
the present investigation and finally, in Sec. IV we give 



a summary of our study. 



II. NONLINEAR GROSS-PITAEVSKII 
EQUATION 

At zero temperature, the time-dependent Bose- 
Einstein condensate wave function ^'(r; r) in two dimen- 
sions at position r and time r may be described by the 
following mean-field nonlinear GP equation jl^, |l^ 



— V2 + V{r)+gN\^ir;T)\'' -flL, 
2m 



OT 



^'(r;T) = 0. 



(2.1) 



Here m is the mass and N the number of atoms in the 
condensate, g the strength of interatomic interaction. 
The trapping potential in two dimensions may be writ- 
ten as V{y) = imu;^(i^ + y^) where lo is the angular 

frequency. The rotational term —ClLz = ihfl{xdy — ydx) 
where Cl is the rotational frequency along the z direc- 
tion with Lz the z component of angular momentum. 
The normalization condition of the wave function is 
J c?r|^(r;r)p = 1. The GP equation (2.1) can accom- 



modate quantized vortex states with rotational motion 
of the condensate around the z axis. 

It is convenient to use the dimensionless variables de- 
fined by a; = x/l, y = y/l, t = tuj/2, I = ^/hJJrmJ), 
il = 2il/LU, and Lp{x,y;t) = l'i'{r;t). In these units Eq. 
(p.l[) becomes 



dx"^ 



dy^ 



M\(f{x,y;t)\ 



dt 



¥>(x,y;i) -0, (2.2) 



where M = 2gNm/li' . The normahzation condition of 
the wave function is 



oo ^ — oo 



dx / dy\ip{x,y;t)\ = 1. 



(2.3) 



The explicit inclusion of the rotational energy in the 
GP hamiltonian above simulates the time evolution of 
a rotational BEC which permits multiple-vortex forma- 
tion. The initial input solution to the GP equation is 
a circularly symmetric state obtained with il = 0. The 
method of solution with ^ is time iteration using an 
appropriate algorithm. If one iterates the GP equation 
above in time using the full angular momentum term the 
circular symmetry is maintained always and a solution 
with multiple vortices is not generated. It is essential 
that the circular symmetry of the solution is broken at 
some stage before multiple vortices are generated during 
time iteration. 

The circular symmetry of the solution can be broken 
in the process of time evolution of the GP equation above 
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by introducing p4| a phase in the initial wave function 
for = of the form 



greater than a critical value. 



S = ^arctan[(y - yi)/{x - Xi)], 



(2.4) 



III. NUMERICAL RESULTS 



where Xi and yi are a set of parameters such that the 
positions (xi,yi) fall inside the boundary of the BEC. 
This initial phase preserves the norm of the wave function 
but accelerates the formation of vortex nucleation after 
time evolution of the GP equation. However, there is 
no direct relation between the parameters {xi,yi) and 
the positions or the numbers of the vortices in the final 
BEC. Another way of facilitating vortex formation is to 
include a small absorptive term in the GP equation as 
has been advocated in Ref. and employed in Ref. 

We solve the GP equation (2.2) with = using 
a time-iteration method elaborated in Refs. |3^ 
starting from the analytic two-dimensional harmonic- 
oscillator solution for TV = 0. The full GP Hamiltonian is 
conveniently broken into three parts H^, Hy and ifres ~ 
the first containing the x-dependent derivative terms, the 
second containing the y-dependent derivative terms with 
the harmonic-oscillator and nonlinear interaction terms 
contributing to the third residual part. The GP equa- 
tions along the x and y directions are discretized on a 
two-dimensional set of grid points x Ny using the 
Crank-Nicholson discretization method |^ . The second 
(first) derivatives are discretized using a three-point (two- 
point) finite difference rule. Effectively, each time iter- 
ation of the GP equation is broken up into three parts 
— using Hj;, Hy and flres- The first two parts lead to 
two sets of tridiagonal equations along x and y direc- 
tions which are solved alternately by the Gaussian elimi- 
nation method ]3l| ] . The third part does not contain any 
derivative and is solved essentially exactly. For a small 
time step A the error involved in this break-up procedure 
along X and y directions is quadratic in A and hence can 
be neglected. Using the above procedure the solution of 
the GP equation corresponding to = and AA / is 
obtained. This solution serves as the initial input for the 
simulation of multiple vortex generation using the time 
iteration of the GP equation. 

The initial solution with O = is circularly symmetric. 
For the generation of multiple vortex the circular sym- 
metry is broken by introducing the phase S of Eq. (2.4) 
to the solution with = and iterating the GP equation 
in time with full fi. We also included a small dissipative 
term in the GP equation. As the formation and decay 
of vortices is intrinsically dissipative l^jjt is more realis- 
tic to include a small dissipative term [g6[ However, 
from a mathematical point of view, vortices and vortex 
lattice can be generated in numerical simulation without 
including any dissipative term |24|. The time evolution 
of the GP equation is performed as described above using 
three independent steps. The system passes through dy- 
namically unstable deformed configurations and during 
time iteration the multiple vortex centers appear for 



In the present numerical study on vortices we use an 
equilibrated condensate with M = 500 trapped in the sta- 
tionary potential as the initial state. Figure 1 shows the 
dynamics of the condensate as it starts to rotate suddenly 
with an angular frequency = 1.5. In the simulation we 
use as input the wave function with phase S given by 
Eq. (2.4). In the beginning the system passes through 
deformed configurations with unstable oscillating bound- 
aries. The vortices are generated in the peripheral region 
via the oscillation of the boundary. More and more vor- 
tices are generated which move towards the center and 
which arrange themselves in a regular form in the con- 
densate. For a sufficiently large number of vortices they 
are arranged as in a lattice and maintain that pattern 
during rotation. When such a regular equilibrated final 
pattern is formed the condensate attains the symmetry 
of the problem — circular in the present case. 





FIG. 1: Formation of the vortices for Af = 500 and fl — 1.5 
at different times; size of each square is 10 x 10. 

In actual experiment the size of the condensate reduces 
with time due to inelastic atomic collisions [||. However, 
throughout the present study we renormalize the number 
of particles during time iteration so that the size is con- 
served. This should not affect the generation of vortices. 

For n = 1.5 and M = 500 one has a hexagonal lattice 
structure with seven vortices at large times. In recent 
theoretical studies |2^, |2^ some general features of 
vortex formation have been established which are con- 
firmed in the present study. For a fixed large TV there 
is a critical 51 for the appearance of a certain number of 
vortices ^, |8|. The critical value of increases as the 
number of vortices increases. For a fixed f2 the number 
of vortices increases as TV increases. It is interesting to 
recall that for TV = only states with a single vortex at 
the center can be formed. 

An interesting set of experiments can be performed 
where on a stable rotating BEC with vortices the strength 
of atomic interaction or the harmonic oscillator trapping 
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potential is suddenly changed. Once either of these pa- 
rameters are changed suddenly on a vortex lattice with 
hexagonal structure, the system enters again a dynami- 
cally unstable asymmetric configuration with highly os- 
cillating deformed boundary. Through the oscillating de- 
formed boundary new vortices can enter the condensate 
or some of the existing vortices can get out. Eventu- 
ally, a new equilibrated BEG with vortices is formed in 
the condensate with circular symmetry. However, as the 
condensate oscillates after the application of the impul- 
sive force it is more difficult to obtain perfect symmetry 
numerically in the final state. Small imperfections re- 
main near the boundary and in the arrangement of the 
vortices in the final state. 
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FIG. 2: Dynamics of the condensate of Fig. 1 at different 
times t upon sudden increase of Af from 500 to 1200 at t = 0; 
size of each square is 10 x 10. 




FIG. 3: Dynamics of the condensate of Fig. 1 at different 

times t upon sudden decrease of A/" from 500 to 100 at f = 0; 
size of each square is 10 x 10. 

In the simulation presented in Fig. 2 the nonlinear 
parameter M on the hexagonal vortex lattice of Fig. 1 
is suddenly changed to 1200 from 500 at t = and the 
evolution of the resultant condensate studied. This cor- 
responds to an increase of the strength of interatomic 



interaction by a factor 2.4. The size of the condensate 
increases and it accommodates more vortices. The sys- 
tem seems to keep on oscillating with (12 ± 1) vortices. 
This oscillation seems to be a natural consequence of the 
sudden change of the strength of interaction. Next we 
consider in Fig. 3 the sudden change in the nonlinear 
parameter J\f from 500 to 100 at t = on the final con- 
figuration of Fig. 1. With the reduction of the strength of 
interaction by a factor 0.2, the condensate shrinks in size 
by passing through dynamically unstable and asymmet- 
ric configurations and starts to oscillate. The condensate 
tries to reach a equilibrated stage by getting rid of vor- 
tices. However, at large times the system continues to 
oscillate with (6 ± 1) vortices. 




FIG. 4: Dynamics of the condensate of Fig. 1 at different 

times upon sudden increase of the harmonic oscillator trap- 
ping potential by a factor of 1.5 at i = 0; size of each square 
is 10 X 10. 




FIG. 5: Dynamics of the condensate of Fig. 1 at different 
times upon sudden decrease of the harmonic oscillator trap- 
ping potential by a factor of 0.75 at t = 0; size of each square 
is 10 X 10. 

Next we consider a sudden change in the harmonic os- 
cillator potential. First, the harmonic oscillator potential 
is increased by a factor of 1.5 on the rotating BEG of Fig. 
1 at t = 0. The time evolution of the resultant conden- 
sate is shown in Fig. 4. With the sudden increase in the 
trapping potential the condensate starts to shrink and 
oscillate. It passes through dynamically unstable asym- 
metric configurations and reaches an equilibrium stage 
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FIG. 6: Decay of the condensate of Fig. 1 at different times 
after the rotation is suddenly stopped at f = 0; size of each 
square is 10 x 10. 




FIG. 7: Free expansion of the condensate Fig. 1 at different 
times after the removal of the trap at t = 0; size of each square 
is 20 X 20. 



by getting rid of a couple of vortices. Eventually, a con- 
densate emerges with four vortices. The system executes 
small oscillations at large times with a fixed number of 
vortices. In the next simulation the harmonic oscillator 
potential is reduced by a factor of 0.75 on the configura- 
tion of Fig. 1 at t = 0. Because of the sudden decrease in 
the trapping potential the condensate starts to swell and 
oscillate. More vortices are generated near the bound- 
ary which move inwards. Eventually, one has a larger 
condensate with 11 vortices at large times. 

We study the decay of the vortices of Fig. 1 once the 
rotation of the BEC is suddenly stopped at t — 0. The 
dynamics of the decay of vortices has been studied exper- 
imentally in Refs. [^P. The numerical simulation of the 
decay of vortices is shown in Fig. 6. The vortices do not 
decay immediately after stopping the rotation. Rather 
they stay for a long time after the rotation stops. As 
time passes vortices get out of the condensate until one 
with a single vortex results. This BEC with a single vor- 
tex survives for a long time. If we consider lo = 2tt x 219 
s~^ as in the experiment of Madison et al one unit of 
present dimensionless time is 0.000726 s. The system gets 
rid of six vortices in about 600 units of time or in 0.4 s. 
However, the state with one vortex stays up to t w 3000, 
which gives a lifetime of about 2400 units of time or about 
1.6 s. Madison et al measured the half life of a single 
vortex state and found it to be about 1 s, which is of the 
same order of magnitude as in the present simulation. 
However, a quantitative comparison between the two is 
not to the point because of approximate nature of the 
present simulation compared to experiment. The present 
calculation used a quasi two-dimensional model opposed 
to a full three-dimensional one. Also, the nonlinearity 
of the present model is not the same as in the experi- 
ment. However, in a rotating BEC with many vortices, 
most of the vortices decay within one second whereas a 
small number of vortices survive for a long time. This 
was also experimentally confirmed by Abo-Shaeer et al 
1^ who used u; « 27r x 85 s^^. 

Finally, we consider the free expansion of the vortices 
in Fig. 1 when the harmonic oscillator trap is suddenly 



removed. The dynamics of free expansion is shown in 
Fig. 7. It is interesting that the vortex pattern is not 
destroyed immediately after the removal of the trap. On 
the other hand, upon expansion the radius of the vortex 
core is found to increase more rapidly than the radius of 
the condensate. Due to the very small radius of the vor- 
tex core it is difficult to observe the vortices and count 
them properly in an experiment |^. However, because of 
the increase in the radius of the vortex core after a free 
expansion the detection and counting of vortices will be 
easier from a picture of the condensate after expansion 
Q. In Fig. 7 the vortices are more prominently visible 
at t = 1.6 than at t = 0.2. Using the radial trapping 
frequency = 27r x 85 s^^ as in the experiment of Abo- 
Shaeer et al ||], i = 1.6 of Fig. 7 corresponds to about 3 
ms. In the actual experiment |9j after 42.5 ms of expan- 
sion the vortex cores were magnified 20 times their size 
in the trap. The present simulation has again led to the 
right order of magnitude. 



IV. SUMMARY 

We have investigated the effect of an impulsive force 
on a rotating BEC with vortices using the mean-field CP 
equation. We have considered the effect of suddenly in- 
creasing or decreasing the strength of atomic interaction 
and the harmonic oscillator trapping potential. Upon 
application of such an impulsive force the condensate 
passes through highly unstable distorted configurations 
with oscillating boundary. New vortices can enter the 
condensate via the oscillating boundary. Also the con- 
densate can get rid of some of the vortices in a simi- 
lar fashion. Eventually, a (more) symmetric condensate 
with a different number of vortices result. When the 
atomic interaction is suddenly increased (decreased) the 
final condensate grows (shrinks) in size and accommo- 
dates more (less) vortices. When the trapping potential 
is suddenly increased (decreased) the final condensate 
shrinks (grows) in size and accommodates less (more) 
vortices. Experiments can be performed when similar 
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impulsive force in applied on a rotating BEC with vor- 
tices and compared with the mean-field prediction. 

We also study the decay of vortices in a rotating BEC 
when the rotation is suddenly stopped. Most of the vor- 
tices decay quickly (in a fraction of a second) , but a few 
vortices survive during a relatively long time before they 
decay. In a typical experimental situation ^ the single 
vortex at the end survives for a unusually long interval 
of time (a couple of seconds) . 

In another simulation of a free expansion of a rotating 
BEC we find that the vortices are not immediately lost 
during expansion. Rather, upon expansion the radius 
of the vortex core increases at a faster rate than that of 
the condensate. Hence the vortices are more prominently 
visible in a picture of the condensate during the free ex- 
pansion of the condensate. Experimentally, it is difhcult 
to observe and count the vortices correctly as the vortex 



cores have usually small radii ||9| . The present simulation 
suggests that the free expansion could aid in the observa- 
tion and counting of vortices in a rotating BEC, specially, 
when there is a large number of vortices. This was noted 
in the experiment of Abo-Shaeer et al ^ . In a typical ex- 
perimental situation, on free expansion the vortices may 
be observed during few tens of milliseconds. 
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